﻿using System;
using System.IO;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.Threading;

namespace SpectrumToWaveLength
{
    /// <summary>
    /// yyyyMMddHHmmss.dat 为光谱文件
    /// yyyyMMddHHmmss.wl 为波长文件
    /// yyyyMMddHHmmss.tmplx 为临时波长文件
    /// </summary>
    public class AnalysingFile
    {
        #region  把光谱数据转换为波长数据

        private const int BUFFERSIZE = 32 * 1024 * 1024;

        public void ProcessFile(object o)
        {
            FileInfo CurFile = o as FileInfo;
            try
            {
                //文件不存在，返回
                if (CurFile == null) return;

                string strName = CurFile.FullName.ToLower().Replace(".dat", ".tmplx");
                //文件已分析，返回
                if (File.Exists(strName)) return;

                //读取原始文件
                FileStream _fsSrc = File.OpenRead(CurFile.FullName);
                BinaryReader _brSrc = new BinaryReader(_fsSrc);
                _fsSrc.Seek(0, SeekOrigin.End);
                byte[] bytes = null;
                //块数据生成的波长文件
                Queue<double[]> qTmpWavLens = new Queue<double[]>(BUFFERSIZE);

                //生成定长波长文件
                FileStream _fsDest = File.OpenWrite(strName);
                _fsDest.SetLength(_fsSrc.Length / 4L);
                _fsDest.Close();

                //从后往前读取数据块
                while (_fsSrc.Position > 0)
                {
                    if (_fsSrc.Position > BUFFERSIZE)
                    {
                        _fsSrc.Seek(-BUFFERSIZE, SeekOrigin.Current);
                        bytes = _brSrc.ReadBytes(BUFFERSIZE);
                        _fsSrc.Seek(-BUFFERSIZE, SeekOrigin.Current);
                    }
                    else
                    {
                        int lentmp = (int)_fsSrc.Position;
                        _fsSrc.Seek(-lentmp, SeekOrigin.Current);
                        bytes = _brSrc.ReadBytes(lentmp);
                        _fsSrc.Seek(-lentmp, SeekOrigin.Current);
                    }
                    AnalysizingBlock(bytes, ref qTmpWavLens);
                    if (qTmpWavLens.Count > 0)
                    {
                        byte[] wavBuffer = new byte[128 * qTmpWavLens.Count];
                        double[][] dd = qTmpWavLens.ToArray();
                        qTmpWavLens.Clear();
                        for (int m = 0; m < dd.Length; m++)
                        {
                            for (int n = 0; n < 8; n++)
                            {
                                Array.Copy(BitConverter.GetBytes(dd[m][n]), 0, wavBuffer, m * 128 + n * 8, 8);
                            }
                        }
                        long lCurStartIndex = _fsSrc.Position / 4L;
                        FileStream fsDest = File.OpenWrite(strName);
                        fsDest.Seek(lCurStartIndex, SeekOrigin.Begin);
                        fsDest.Write(wavBuffer, 0, wavBuffer.Length);
                        fsDest.Close();
                    }

                }
            }
            catch (Exception ex) { System.Diagnostics.Debug.WriteLine(ex.StackTrace + ":" + ex.Message); }
        }
        
        /// <summary>
        /// 处理记录文件块
        /// </summary>
        /// <param name="BufferSrc">记录文件块</param>
        /// <param name="queWavArray">波长记录</param>
        private void AnalysizingBlock(byte[] BufferSrc, ref Queue<double[]> qCurWaveLenths)
        {
            qCurWaveLenths = new Queue<double[]>();
            try
            {
                double[] _wavlen = new double[16];
                double[] _tmpwavlen = new double[16];
                double[] _tmpData = new double[256];
                int iCurRecordsCount = BufferSrc.Length / 512;
                for (int i = 0; i < iCurRecordsCount; i++)
                {
                    Array.Clear(_tmpData, 0, 256);
                    Array.Clear(_wavlen, 0, _wavlen.Length);
                    for (int n = 0; n < 256; n++)
                    {
                        _tmpData[n] = double.Parse((BitConverter.ToUInt16(BufferSrc, 512 * i + 2 * n)).ToString());
                    }
                    double dAvrg = _tmpData.Average() + 2000;
                    for (int j = 0; j < _tmpData.Length; j++)
                    {
                        if (_tmpData[j] <= dAvrg)
                        {
                            _tmpData[j] = 0;
                        }
                    }
                    GetWavesLength(_tmpData, ref _wavlen);
                    Array.Clear(_tmpwavlen, 0, 16);
                    if (_wavlen != null)
                    {
                        if (_wavlen.Length < 17)
                        {
                            Array.Copy(_wavlen, 0, _tmpwavlen, 0, _wavlen.Length);
                        }
                        else
                        {
                            Array.Copy(_wavlen, 0, _tmpwavlen, 0, 16);
                        }
                    }
                    qCurWaveLenths.Enqueue((double[])_tmpwavlen.Clone());
                }
            }
            catch (Exception ex) { System.Diagnostics.Debug.WriteLine(ex.StackTrace + ":" + ex.Message); }
        }
        
        #endregion

        #region 获取波长值
        private static double A = 1570.398804;
        private static double B1 = -0.1617357432842254600000000;
        private static double B2 = -0.0000663390237605199220000;
        private static double B3 = -0.0000000354414062542218740;
        private static double B4 = 0.0000000001906033247900396;
        private static double B5 = -0.0000000000003401445791695;
        private static double[] PreWav = new double[32];

        private static void GetWavesLength(double[] Spectra, ref double[] wavelength)
        {
            int fbgNo = 0;
            GetWavesLength(Spectra, ref wavelength, ref fbgNo);
        }
        /// <summary>
        /// 
        /// </summary>
        /// <param name="Spectra"></param>
        /// <param name="wavelength"></param>
        /// <param name="fbgNo"></param>
        private static void GetWavesLength(double[] Spectra, ref double[] wavelength, ref int fbgNo)
        {
            List<double> lstwav = new List<double>();
            fbgNo = 0;
            int iCount = 0;
            double[] _x = new double[256], _y = new double[256];
            for (int i = 0; i < Spectra.Length; i++)
            {
                if (Spectra[i] > 0)
                {
                    _x[iCount] = i;
                    //_y[iCount++] = Math.Log(Spectra[i]);
                    _y[iCount++] = Math.Pow(Spectra[i], 2);
                }
                else if (iCount > 0)
                {
                    if (iCount > 2)
                    {
                        double[] _tx = new double[iCount];
                        double[] _ty = new double[iCount];
                        while (iCount-- > 0)
                        {
                            _tx[iCount] = _x[iCount];
                            _ty[iCount] = _y[iCount];
                        }
                        //PreWav[fbgNo] = _wavlen(_tx, _ty);
                        PreWav[fbgNo] = (PreWav[fbgNo] > 0) ? Avrdouble(PreWav[fbgNo], _wavlen(_tx, _ty)) : _wavlen(_tx, _ty);
                        if (PreWav[fbgNo] < 1525 || PreWav[fbgNo] > 1570)
                        {
                            continue;
                        }
                        lstwav.Add(Math.Round(PreWav[fbgNo++], 5));
                    }
                    _x = new double[256];
                    _y = new double[256];
                    iCount = 0;
                }
            }
            lstwav.Sort();
            wavelength = lstwav.ToArray();
        }
        private static double Avrdouble(double _dpre, double _d)
        {
            double _t = Math.Abs(Math.Floor((_d - _dpre) * 1000d));
            if (_t < 180) return _d;
            if (_t < 250) return _dpre;
            return 0.25d * _d + 0.75d * _dpre;
        }

        /// <summary>
        /// 通过高斯拟合（最小二乘法求参数）求出波长值
        /// </summary>
        /// <param name="_X">因子</param>
        /// <param name="_Y">因子</param>
        /// <returns>波长值</returns>
        private static double _wavlen(double[] _X, double[] _Y)
        {
            if (_X.Length != _Y.Length || _Y.Length < 3)
            {
                return 0;
            }
            try
            {
                double s00 = 0, s10 = 0, s20 = 0, s30 = 0, s40 = 0, s01 = 0, s11 = 0, s21 = 0;
                for (int i = 0; i < _X.Length; i++)
                {
                    s00 += 1;//sigma x^0*y^0;
                    s10 += _X[i];//sigma x;   
                    s20 += Math.Pow(_X[i], 2);//sigma x^2;
                    s30 += Math.Pow(_X[i], 3);//sigma x^3;
                    s40 += Math.Pow(_X[i], 4);//sigma x^4
                    s01 += _Y[i]; //sigma y;
                    s11 += _X[i] * _Y[i]; //sigma x*y;
                    s21 += _X[i] * _X[i] * _Y[i];//sigma x^2*y;
                }
                double tmp = 0.5d * (s21 * (s20 * s10 - s30 * s00) + s11 * (s40 * s00 - s20 * s20) + s01 * (s30 * s20 - s40 * s10))
                    /
                        (s21 * (s10 * s10 - s20 * s00) + s11 * (s30 * s00 - s20 * s10) + s01 * (s20 * s20 - s30 * s10));

                return (A + B1 * tmp + B2 * Math.Pow(tmp, 2) + B3 * Math.Pow(tmp, 3) + B4 * Math.Pow(tmp, 4) + B5 * Math.Pow(tmp, 5));
            }
            catch (Exception ex) { System.Diagnostics.Debug.WriteLine(ex.StackTrace + ":" + ex.Message); }
            return 0;
        }
        #endregion

    }
}

